﻿using ESRI.ArcGIS.Carto;
using ESRI.ArcGIS.Geodatabase;
using ESRI.ArcGIS.Geometry;
using System;
using System.Collections.Generic;
using System.IO;
using System.Linq;
using System.Text;
using System.Windows.Forms;

namespace TrafficAccessibilityIndex
{
    struct NetPointInf
    {
        public int EdgeID;
        public IPoint edgePoint;
        public double percent;
        public double distance;
    }
    class ParkGravityModel : IGravityModel
    {
        public ILayer m_Layer;
        public ILayer m_PoiLayer;//
        public List<double> m_LayerAccessability;//layer图层各个要素的可达性，key为要素ID
        public List<List<double>> m_LayerNetDistance;//记录layer的每个要素到poiLayer每个要素的网路距离
        public ClsPathFinder m_PathFinder;//网路距离计算类
        public ILayer m_JounctionLayer;//路网节点图层
        public ToolStripStatusLabel statusBar;//工具栏
        /// <summary>
        /// 初始化
        /// </summary>
        /// <param name="layer1">网格图层</param>
        /// <param name="layer2">POI图层</param>
        /// <param name="layer3">路网节点图层</param>
        /// <param name="pathFinder"></param>
        /// <param name="bar">状态栏</param>
        public ParkGravityModel(ILayer layer1, ILayer layer2, ILayer layer3, ClsPathFinder pathFinder, ToolStripStatusLabel bar)
        {
            m_Layer = layer1;
            m_PoiLayer = layer2;
            m_PathFinder = pathFinder;
            m_LayerAccessability = new List<double>();
            m_LayerNetDistance = new List<List<double>>();
            m_JounctionLayer = layer3;
            statusBar = bar;
        }


        public void ReadLayerNetDistance(string path)
        {
            this.statusBar.Text = "正在读取路网距离。";
            Application.DoEvents();
            StreamReader sw = new StreamReader(path);
            string str = sw.ReadLine();
            int i = 0;
            while (str != null)
            {
                if (i % 100 == 0)
                    GC.Collect();
                List<double> values = new List<double>();
                string[] strs = str.Split('_');
                foreach (string item in strs)
                {
                    if (item != "")
                        values.Add(Convert.ToDouble(item));
                }
                m_LayerNetDistance.Add(values);
                Array.Clear(strs, 0, strs.Length);
                str = sw.ReadLine();
                i++;
            }
            sw.Close();
            GC.Collect();
        }

        public void SaveLayerNetDistance(string path)
        {
            this.statusBar.Text = "正在保存路网距离。";
            Application.DoEvents();
            StreamWriter sw = new StreamWriter(path);
            foreach (List<double> lists in m_LayerNetDistance)
            {
                StringBuilder sb = new StringBuilder();
                foreach (double value in lists)
                {
                    sb.Append(value + "_");
                }
                sw.WriteLine(sb);
            }
            sw.Close();
            GC.Collect();
        }

        public void CalcLayerNetDistance(ref StreamWriter sw)
        {
            string filename = "";
            OpenFileDialog ofd = new OpenFileDialog();
            ofd.Title = "读取网络距离文件";
            ofd.Filter = "*.lnd|*.lnd";
            if (ofd.ShowDialog() == DialogResult.OK)
            {
                filename = ofd.FileName;
                ReadLayerNetDistance(filename);
                sw = new StreamWriter(filename, true);
            }
            else
            {
                SaveFileDialog sfd = new SaveFileDialog();
                sfd.Title = "保存网络距离";
                sfd.Filter = "*.lnd|*.lnd";
                if (sfd.ShowDialog() == DialogResult.OK)
                {
                    filename = sfd.FileName;
                    sw = new StreamWriter(filename);
                }
            }
            try
            {
                List<NetPointInf> p_FeatureNPI = new List<NetPointInf>();
                List<NetPointInf> p_PoiNPI = new List<NetPointInf>();
                IFeatureLayer p_FeatureLayer = m_Layer as IFeatureLayer;
                IFeatureLayer p_PoiLayer = m_PoiLayer as IFeatureLayer;
                foreach (var tmp in m_LayerNetDistance)
                {
                    p_FeatureNPI.Add(new NetPointInf());
                }
                int p_PeatureCount = p_FeatureLayer.FeatureClass.FeatureCount(null);
                int p_PoiCount = p_PoiLayer.FeatureClass.FeatureCount(null);
                string tmpStr1 = "共" + p_PeatureCount + "个格网" + p_PoiCount + "个POI对象";
                IFeatureCursor p_FeatureCursor = p_FeatureLayer.FeatureClass.Search(null, false);

                IFeature p_Feature = p_FeatureCursor.NextFeature();
                int i = 1;
                while (p_Feature != null)
                {
                    if (m_LayerNetDistance.Count >= i)
                    {
                        p_Feature = p_FeatureCursor.NextFeature();
                        i++;
                        continue;
                    }

                    StringBuilder sb = new StringBuilder();
                    List<double> p_FeatureNetDistance = new List<double>();
                    m_LayerNetDistance.Add(p_FeatureNetDistance);
                    IPoint p_FeaturePt = GetFeatureCenterPoint(p_Feature);
                    IFeatureCursor p_PoiCursor = p_PoiLayer.Search(null, false);
                    IFeature p_Poi = p_PoiCursor.NextFeature();
                    string tmpStr2 = tmpStr1 + ";正在计算第" + i + "个格网距离第";
                    int j = 1;
                    while (p_Poi != null)
                    {
                        if (((i - 1) * p_PeatureCount + j) % 20 == 0)
                        {
                            statusBar.Text = tmpStr2 + j + "个公园的距离";
                            GC.Collect();
                            Application.DoEvents();
                        }
                        IPoint p_PoiPt = GetFeatureCenterPoint(p_Poi);
                        double netDis = CalcNetDistance(p_FeaturePt, p_PoiPt, i - 1, j - 1, p_FeatureNPI, p_PoiNPI);
                        if (netDis == 0)
                            netDis = 0;
                        p_FeatureNetDistance.Add(netDis);
                        sb.Append(netDis + "_");
                        p_Poi = p_PoiCursor.NextFeature();
                        j++;
                    }
                    System.Runtime.InteropServices.Marshal.ReleaseComObject(p_PoiCursor);
                    sw.WriteLine(sb);
                    sb.Clear();
                    p_Feature = p_FeatureCursor.NextFeature();
                    i++;
                }
                p_FeatureNPI.Clear();
                p_PoiNPI.Clear();
                sw.Close();
                System.Runtime.InteropServices.Marshal.ReleaseComObject(p_FeatureCursor);
            }
            catch (Exception ex)
            {
                sw.Close();
            }

        }

        public double CalcNetDistance(ESRI.ArcGIS.Geometry.IPoint ptCell, ESRI.ArcGIS.Geometry.IPoint ptPoi, int i, int j, List<NetPointInf> p_FeatureNPI, List<NetPointInf> p_PoiNPI)
        {
            MultipointClass ptsForFindPath = new MultipointClass();
            m_PathFinder.StopPoints = ptsForFindPath;
            object o = Type.Missing;
            if (p_FeatureNPI.Count - 1 < i)
            {
                int intEdgeID = int.MinValue;
                double dblEdgePercent = double.MinValue;
                IPoint ipFoundEdgePoint = null;
                m_PathFinder.m_ipPointToEID.GetNearestEdge(ptCell, out intEdgeID, out ipFoundEdgePoint, out dblEdgePercent);
                if (ipFoundEdgePoint == null)
                {
                    ipFoundEdgePoint = GetNearestJounctionPoint(ptCell, m_JounctionLayer);
                    m_PathFinder.m_ipPointToEID.GetNearestEdge(ipFoundEdgePoint, out intEdgeID, out ipFoundEdgePoint, out dblEdgePercent);
                }
                NetPointInf npi = new NetPointInf();
                npi.EdgeID = intEdgeID;
                npi.edgePoint = ipFoundEdgePoint;
                npi.percent = dblEdgePercent;
                npi.distance = Math.Sqrt((ptCell.X - ipFoundEdgePoint.X) * (ptCell.X - ipFoundEdgePoint.X) +
                     (ptCell.Y - ipFoundEdgePoint.Y) * (ptCell.Y - ipFoundEdgePoint.Y));
                p_FeatureNPI.Add(npi);
            }
            ptsForFindPath.AddPoint(p_FeatureNPI[i].edgePoint, ref o, ref o);



            if (p_PoiNPI.Count - 1 < j)
            {
                int intEdgeID = int.MinValue;
                double dblEdgePercent = double.MinValue;
                IPoint ipFoundEdgePoint = null;
                m_PathFinder.m_ipPointToEID.GetNearestEdge(ptPoi, out intEdgeID, out ipFoundEdgePoint, out dblEdgePercent);
                if (ipFoundEdgePoint == null)
                {
                    ipFoundEdgePoint = GetNearestJounctionPoint(ptPoi, m_JounctionLayer);
                    m_PathFinder.m_ipPointToEID.GetNearestEdge(ipFoundEdgePoint, out intEdgeID, out ipFoundEdgePoint, out dblEdgePercent);
                }
                NetPointInf npi = new NetPointInf();
                npi.EdgeID = intEdgeID;
                npi.edgePoint = ipFoundEdgePoint;
                npi.percent = dblEdgePercent;
                npi.distance = Math.Sqrt((ptPoi.X - ipFoundEdgePoint.X) * (ptPoi.X - ipFoundEdgePoint.X) +
                     (ptPoi.Y - ipFoundEdgePoint.Y) * (ptPoi.Y - ipFoundEdgePoint.Y));
                p_PoiNPI.Add(npi);
            }
            ptsForFindPath.AddPoint(p_PoiNPI[j].edgePoint, ref o, ref o);

            try
            {
                int[] intEdgeIDS = new int[2];
                intEdgeIDS[0] = p_FeatureNPI[i].EdgeID;
                intEdgeIDS[1] = p_PoiNPI[j].EdgeID;
                double[] EdgePercent = new double[2];
                EdgePercent[0] = p_FeatureNPI[i].percent;
                EdgePercent[1] = p_PoiNPI[j].percent;
                double tmpDis = p_FeatureNPI[i].distance + p_PoiNPI[j].distance;
                int id = m_PathFinder.SolvePath("Length", intEdgeIDS, EdgePercent);

                IPoint ipFoundEdgePoint = p_FeatureNPI[i].edgePoint;
                IPoint ipFoundEdgePoint1 = p_PoiNPI[j].edgePoint;
                double cost = m_PathFinder.PathCost + tmpDis;
                double distance = Math.Sqrt((ipFoundEdgePoint.X - ipFoundEdgePoint1.X) * (ipFoundEdgePoint.X - ipFoundEdgePoint1.X) +
                                             (ipFoundEdgePoint.Y - ipFoundEdgePoint1.Y) * (ipFoundEdgePoint.Y - ipFoundEdgePoint1.Y)) + tmpDis;
                if (cost < distance)
                {

                    cost = distance;
                }

                ptsForFindPath.SetEmpty();
                return cost;
            }
            catch (Exception ex)
            {
                return 0;
            }

        }

        public void CalcLayerAccessability()
        {
            statusBar.Text = "正在计算公园服务面积。";
            Application.DoEvents();
            List<double> Dvalues = CalcDvalues("合计");

            IFeatureLayer p_FeatureLayer = m_Layer as IFeatureLayer;
            string fieldName = "";
            if (m_PoiLayer.Name.Length > 10)
                fieldName = m_PoiLayer.Name.Substring(0, 10);
            else
                fieldName = m_PoiLayer.Name;
            int IDIndex1 = p_FeatureLayer.FeatureClass.FindField(fieldName);
            IFeatureLayer p_PoiLayer = m_PoiLayer as IFeatureLayer;
            IFeatureCursor p_FeatureCursor = p_FeatureLayer.Search(null, false);
            int p_PeatureCount = p_FeatureLayer.FeatureClass.FeatureCount(null);
            IFeature p_Feature = p_FeatureCursor.NextFeature();
            int i = 0;
            while (p_Feature != null)
            {
                if ((i + 1) % 10 == 0)
                {
                    statusBar.Text = "共" + p_PeatureCount + "个网格，正在计算第" + (i + 1) + "个网格";
                    GC.Collect();
                    Application.DoEvents();
                }
                IFeatureCursor p_PoiCursor = p_PoiLayer.Search(null, false);
                int j = 0;
                IFeature p_Poi = p_PoiCursor.NextFeature();
                double accessValue = 0;//可达性
                while (p_Poi != null)
                {
                    IArea area = p_Poi.Shape as IArea;

                    double netDis = m_LayerNetDistance[i][j];
                    accessValue += area.Area / (Dvalues[j] * netDis);//β=1
                    //accessValue += area.Area / (Dvalues[j] * netDis * netDis);//β=2
                    p_Poi = p_PoiCursor.NextFeature();
                    j++;
                }
                System.Runtime.InteropServices.Marshal.ReleaseComObject(p_PoiCursor);
                p_Feature.set_Value(IDIndex1, accessValue);
                p_Feature.Store();
                m_LayerAccessability.Add(accessValue);
                p_Feature = p_FeatureCursor.NextFeature();
                i++;
            }
            System.Runtime.InteropServices.Marshal.ReleaseComObject(p_FeatureCursor);
        }

        public List<double> CalcDvalues(string name)
        {
            List<double> Dvalue = new List<double>();
            IFeatureLayer p_FeatureLayer = m_Layer as IFeatureLayer;
            IFeatureLayer p_PoiLayer = m_PoiLayer as IFeatureLayer;

            IFeatureCursor p_PoiCursor = p_PoiLayer.Search(null, false);
            IFeature p_Poi = p_PoiCursor.NextFeature();
            int p_FieldIndex = p_FeatureLayer.FeatureClass.FindField(name);
            int i = 0;

            double value = 0;
            while (p_Poi != null)
            {
                int j = 0;
                IFeatureCursor p_FeatureCursor = p_FeatureLayer.Search(null, false);
                IFeature p_Feature = p_FeatureCursor.NextFeature();
                while (p_Feature != null)
                {
                    int population = (int)p_Feature.get_Value(p_FieldIndex);
                    value += population / m_LayerNetDistance[j][i];//β=1
                    //value += population / (m_LayerNetDistance[j][i]*m_LayerNetDistance[j][i]);//β=2
                    p_Feature = p_FeatureCursor.NextFeature();
                    j++;
                }
                System.Runtime.InteropServices.Marshal.ReleaseComObject(p_FeatureCursor);
                Dvalue.Add(value);
                p_Poi = p_PoiCursor.NextFeature();
                i++;
            }
            System.Runtime.InteropServices.Marshal.ReleaseComObject(p_PoiCursor);
            return Dvalue;
        }

        private IPoint GetFeatureCenterPoint(IFeature f)
        {
            if (f.Shape.GeometryType == esriGeometryType.esriGeometryPoint)
                return f.Shape as IPoint;
            else if (esriGeometryType.esriGeometryPolygon == f.Shape.GeometryType)
                return (f.Shape as IArea).Centroid;
            else
                return null;
        }

        private IPoint GetNearestJounctionPoint(IPoint ptCell, ILayer layer)
        {
            IFeatureLayer feaLayer = layer as IFeatureLayer;

            double baseDistance = 500;
            ISpatialFilter spatialFilter = new SpatialFilterClass();
            spatialFilter.Geometry = GetGeometry(ptCell, baseDistance);
            spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelContains;

            IFeatureCursor cursor = feaLayer.FeatureClass.Search(spatialFilter, false);
            IFeature feature = cursor.NextFeature();
            while (feature == null)
            {
                baseDistance += 500;
                spatialFilter.Geometry = GetGeometry(ptCell, baseDistance);
                spatialFilter.SpatialRel = esriSpatialRelEnum.esriSpatialRelContains;
                System.Runtime.InteropServices.Marshal.ReleaseComObject(cursor);
                cursor = feaLayer.FeatureClass.Search(spatialFilter, false);
                feature = cursor.NextFeature();
            }

            double min = double.MaxValue;
            IPoint result = null;
            while (feature != null)
            {
                IPoint pt = feature.Shape as IPoint;

                double distance = Math.Sqrt((ptCell.X - pt.X) * (ptCell.X - pt.X) +
                                             (ptCell.Y - pt.Y) * (ptCell.Y - pt.Y));
                if (min > distance)
                {
                    min = distance;
                    result = pt;
                }

                feature = cursor.NextFeature();
            }
            System.Runtime.InteropServices.Marshal.ReleaseComObject(cursor);
            return result;
        }
        private IGeometry GetGeometry(IPoint pt, double distance)
        {
            ITopologicalOperator topolOperator = pt as ITopologicalOperator;
            return topolOperator.Buffer(distance);
        }
    }
}
